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A linear network analysis is developed for the system design of split-type free-dis¬ 
placer Stirling refrigerators. The dynamics models are derived to describe the 
input/output relation of each component by use of the governing equations in conjunc¬ 
tion with linearization and approximation. By connecting the equivalent circuits of the 
components together, a linear network consisting of a fluid network and a displacer 
network is obtained. The fluid network accounts for the propagation of gas flow 
induced by the pressure wave; while the displacer network accounts for the displacer 
motion induced by the pressure force exerted on the displacer. Two transfer functions 
are further derived for the displacer motion and the gas pressure of the expansion 
space with respect to the piston motion. The system performance is then evaluated 
from the frequency response of the transfer functions by using the sinusoidal signal 
analysis. The performance prediction of a split-type free-displacer Stirling refrigerator 
using the linear network analysis is shown to be satisfactory. © 1996 Elsevier 
Science Limited 

Keywords: Stirling refrigerator; split type Stirling refrigerator; linear network model 
of Stirling cooler 


Nomenclature 

V 

volume (m 3 ) 



X 

distance (m) 

a 

gain 

X 

displacement (m) 

A 

area (m 2 ) 

T 

time constant (s) 

C 

flow capacitance (m’kg/kJ) 

W 

angular frequency (=2irf) 

C d 

displacer loss coefficient (N s/m) 

P 

density 

Cp 

gas specific heat at constant pressure 
(kJ/kg K) 

e 

piston angle (deg) 

C v 

h 

gas specific heat at constant volume (kJ/kg K) 
convection heat transfer coefficient (W/m 2 K) 

Subscripts 

d 

diameter (m) 

c 

compression chamber; characteristic 

f 

frequency (Hz) 


impedance 

k 

thermal conductivity (W/m K) 

d 

displacer 

K u 

gas spring constant (N/m) 

e 

expansion space 

L 

length (m) 

f 

gas; fluid; flow 

m 

mass flow rate (kg/s) 

g 

gas 

R 

gas constant (kJ/(kgK)) 

h 

hot end; hydraulic 

Re 

Reynolds number, dimensionless 

H 

hot end of the regenerator 

P 

pressure (N/m 2 ) 

i 

inlet 

P 

pressure (N/m 2 ) 

L 

cold end 

Pch 

charge pressure (N/m 2 ) 

m 

mean 

Qnct 

net cooling capacity (W) 

n 

natural frequency 

S 

Laplace variable, complex 

0 

outlet; equilibrium point; amplitude 

T 

temperature (K) 

P 

piston 

1 

time (s) 

r 

regenerator 
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s solid material 

t connecting tube 

w warm space 


Superscripts 

steady-state point; equilibrium point 
perturbed value 
time change rate 


The thermal performance analysis of a Stirling refrigerator 
is very important at the hardware design stage. Researchers 
have attempted to develop an analytical method and a com¬ 
puter-aided tool for the design application. But the achieve¬ 
ment is limited so far. The designer still relies heavily on 
personnel experience or trial and error. For a split-type free- 
displacer Stirling refrigerator, the losses of wall friction and 
clearance seal leakage of the displacer are rarely predicted 
accurately. Hence, the analysis can deviate greatly from the 
test results. 

Nodal analysis' is a typical method developed for free- 
displacer Stirling refrigerator simulation. However, the 
time-variant conservation equations of mass, momentum 
and energy have to be solved at each node of the whole 
refrigerator until the solution converges at each time 
instant. The computation is then repeated for different 
designs or operating conditions in order to obtain a per¬ 
formance map of the refrigerator. The solution procedure 
is so complicated and time consuming that it usually 
requires use of a super-computer. In order to simplify the 
analysis and increase the computational speed, a linear net¬ 
work analysis is developed in the present study. 

A split-type free-displacer Stirling refrigerator comprises 
a compressor and a cold finger connected by a connecting 
tube. The cold finger contains a regenerator/displacer, a gas 
spring, and an expansion space (cold head) as shown in 
Figure 1. The displacer is driven by the gaseous pressure 
wave resulting from the reciprocating motion of the piston. 
The whole refrigerator thus operates at a cyclically steady 
state very close to sinusoidal. The refrigerator can thus be 
treated as a dynamical system with the piston displacement 
as the system input and the displacer displacement and the 
expansion space pressure as the system output 2 . 

From the viewpoint of system dynamics, any component 
of a dynamical system is triggered by an input (physical 


force) and induces an output (dynamic response). A 
dynamics model can be derived to describe the input/output 
relation of each component by use of the governing equa¬ 
tions in conjunction with linearization and approximations. 
A network model can be further derived from the electric 
circuit analogy. An equivalent circuit or system block dia¬ 
gram can then be used to represent the dynamic behaviour 
of each component. Connecting the analogous circuits of 
all the components together will lead to a linear network 
of the refrigerator 3 . For the split-type free-displacer Stirling 
refrigerator, a transfer-function model representing the per¬ 
formance of the refrigerator can be derived, and from this 
the system performance can be easily evaluated. A brief 
summary of the analysis has been published in the previous 
technical note 2 . A detailed description of the linear network 
analysis for the split-type free-displacer Stirling refrigerator 
is reported in the present paper. 


Dynamics model of components 

The split-type free-displacer Stirling refrigerator essentially 
consists of six components: namely, compressor, con¬ 
necting tube, warm space, regenerator, displacer, and 
expansion space. The dynamics model of each component 
can be derived. For simplification, the ideal gas assumption 
is used throughout the derivation. 

Compressor 

The compressor of a Stirling refrigerator consists of a com¬ 
pression space and a piston driven by a crank mechanism. 
Since the piston motion is very fast and the mixing in the 
compression chamber is severe, the gas temperature and 
pressure inside the compression chamber are assumed to 
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Figure 1 Schematic diagram of split-type free-displacer Stirling refrigerator 
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be uniform. The actual gas compression process inside the 
compression chamber lies between isothermal and adia¬ 
batic. Since cooling is provided for the compressor of a 
Stirling refrigerator, the isothermal compression with a con¬ 
stant gas temperature T c is assumed. Assuming no gas leak¬ 
age between the piston and the cylinder, conservation of 
mass to the compression chamber will lead to the dynam¬ 
ics model 


m c (t) = - 


RTr 




dp At) 
dt 


-ApPM 


dX p (r) 


dt 


(I) 


m c is the mass flow rate out of the compressor; p c is the 
gas pressure; V c (t) = V co - ApX p (r);V co is the volume of 
compression chamber for the piston at the equilibrium pos¬ 
ition (X p = 0); X p (t) is the piston displacement measured 
from the midpoint of the piston motion toward the top dead 
end; A p is the cross section area of the piston. 

Applying small perturbation around the equilibrium 
point,_i.e. X p (t) = X p + X p (t); p c (t) = p c + p c (t): m c (t ) = 
m c + ih c (t) - ih c (t) to Equation (1), neglecting higher-order 
terms, taking that X p = 0 is the piston central position, and 
then taking Laplace transforms, we obtain a perturbed 
dynamics model 


variation. From the above assumptions, we obtain the gov¬ 
erning equations of the connecting tube from the conser¬ 
vations of mass and momentum. 


1 dp(x,t) 1 dm (x,t) 
RT t dt A, dx 


1 diii(x,r) 

A, dx + 


dp(x,t) 

dx 


K 

+ — m(x,t) = 0 

A t 


(3) 

(4) 


where A, is the tube cross section area. It can be easily 
shown that Equations (3) and (4) are the wave equations. 

Applying small perturbation around the equilibrium point 
and with ih = 0 for a cyclically-steady operation, i.e. 
p(x,t) - p + p(x,t); A(x,t) = th + fh(x,t) = fh(x,t), to Equa¬ 
tions (3) and (4) and then solving the Laplace transformed 
equations, we obtain the dynamics model of the con¬ 
necting tube 


m to (s) 


cosh(T t L t ) 

- - 1 sinh(T t L,) 

Px ois) 


- Z cl sinhf I' t L,) 

^Ct 

cosh(T t L t ) 


r spoA„ ~ ~ p c (s ) 

^ (S)= wr w = ra:(s) + T j ^ c c (2) 

where C c = V C JRT C is the equivalent capacitance of the 
compression chamber; m p (s) acts as a current source rep¬ 
resenting the available or gross mass flow generated by the 
piston motion. The equivalent circuit is shown in Figure 2. 

Connecting tube 

The connecting tube links the compressor and the cold 
finger. The flow field is assumed to be 1-D. In order for 
linearization, the second-order viscous and the inertia terms 
in the momentum equation are neglected. However, the 
second-order viscous friction is taken into account by a 
modified resistance coefficient k which is calculated using 
piecewise linear approximation 4 . Therefore, the viscous 
resistance is proportional to the mass flowrate with a pro¬ 
portional constant k depending upon the amplitude of the 
mass flow rate. 

The gas in the connecting tube is assumed to undergo 
an isothermal process with a constant temperature 7j. This 
can hold since the connecting tube is usually small in diam¬ 
eter and is thick in wall thickness in order to withstand the 
high gas pressure. Hence, the connecting tube can act as 
an energy storage medium to damp out the gas temperature 


m c 





Figure 2 Equivalent circuit of compression chamber 



X 


WtiO) 

Pti( s ) 


(5) 


where T t = yZ, Ij and Z ct = Jz t /Y t are the propagation con¬ 
stant and the characteristic impedance of the tube, respect¬ 
ively; Z t = R Fi + sL Ft = series impedance and V, = ,vC Fl = 
shunt admittance; F, and Z c , are derived as 

T, = yjsCp t (R Fl + .v/. F i); Z ct = (6) 

vc Ft 

where C Ft , L Fl and are the flow capacitance, flow induct¬ 
ance, and flow resistance per unit tube length, respectively, 
which are defined as 


C Ft — 




K 

A t 


(7) 


where 


K = 0.1556(p t w max d t /p,)-°' 201 ( w m Jd h ), (8) 

w max is the peak velocity of the oscillating flow in the tube; 
d, is the hydraulic diameter of the tube 5 . 

It is worth noting that k is considered to be a constant 
during the modelling, but should be adjusted by numerical 
iteration during the computation to give a correct value for 
the corresponding mass flow u' m . LK . 

Since the dynamics model of the connecting tube belongs 
to a distributed-parameter system, an equivalent circuit con¬ 
sisting of infinite sets of shunt and series impedances can 
be drawn as shown in Figure J. Figure 3 approximates the 
model of Equation (5) as N is large. The shunt and series 
impedances for each segment of the tube satisfy the follow¬ 
ing relation 
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-C 


n+l 


Pto 


■O 


- Z 3~ "• ~ Z n- ’ Z l~ Z " +1 “ 2 ( 

Yj = Y 2 = Y 3 = - =Y n = Y t (^) 

Figure 3 Equivalent circuit of connecting tube 


Zi = Z 3 = ■■■ 


-z = z 

'-I 



Z, 


5 

2 



K 1 = F 2 = Y 3 = '-- = y n = y t 



(9) 


The limiting case N & corresponds to the present model. 
Equation (5). A similar circuit can also be drawn based on 
the series expansion of coshfF,/-,) and sinh(r,L,) in Equ¬ 
ation (5) with respect to r,L t . 

A system block diagram as shown in Figure 4 is used 
here to illustrate the input/output relation of the connecting 
tube from the system dynamics point of view. In practice, 
Figure 3 will be used as the analogous circuit during the 
linear network analysis. 

From Equation (5), we also obtain Jhe relation between 
the output impedance Z to _ [ = p ta (s)/m ta (s)] and the input 
impedance Z ti [ = /7 ti (^)/m ti (s)] of the connecting tube 

1 + ^' tanh[ r,L, | 

Z ti =- z ---Z to (10) 

1 + -- tanh[r,(s)L,] 

Zj C \ 


Warm space 

The gas temperature is assumed uniform in the warm space 
since the gas mixing is severe. Similar to the compression 
chamber, a linearly perturbed model can be derived from 
mass continuity equation 



Figure 4 Block diagram of connecting tube 


~ , , ~ , , . ■sp^dw , , Pw( s ) ,,,, 

m v ,(s) = m wo (s) + - —— X d (s) + —— (11) 

l/.vC w 

The second term in the right-hand side of Equation (11) is 
the mass flow rate induced by the displacer motion, which 
is analogous to electric current source; the third term rep¬ 
resents the flow capacitance effect of the warm space; C w 
= VJ( RT^). The equivalent electric circuit is shown in 
Figure 5. 

Regenerator 

The regenerator is an energy-storage element made from 
wire screens. The derivation of the dynamical model is 
similar to that of the connecting tube. 

For the momentum equation, the inertia term 
(l/Af r )d(m|m|/p)/dx and the second-order viscous term 
pl$(dpA fr ) 2 m\m\ can be neglected. This can hold since the 
Reynolds number in the regenerator is not large. The press¬ 
ure loss due to the second-order viscous friction is taken 
into account by a modified frictional coefficient a calcu¬ 
lated using piecewise linear approximation 4 . Therefore, the 
viscous resistance is assumed to be proportional to the mass 
flow with a proportional constant a which depends upon 
the amplitude m max of the oscillating flow, a is considered 
to be a constant during the modelling, however it should 
be adjusted by numerical iteration during the computation 
in order to give a correct value for the corresponding mass 
flow (m max ). 

Assuming 1-D flow, no axial conduction and constant 
properties, the transient governing equations in terms of 
perturbed variables m(x,t) =_m r (x) + >h(x,t).p(x,t) = P T {x) + 
p(x,t), and T(x,t) = T r (x) + T(x,t) and using in X*) = 0 for 
cyclically-steady operation, are derived from the conser¬ 
vation of mass and momentum 


rh W i( s ) rh W o(s) 



Figure 5 Equivalent circuit of warm space 
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Continuity equation of gas 


T( , dP(x ’H - dT(x,t) 
T t (x) —-- P r (x) 


at 

RT?(x) dm(x,t) 
A fr dX 


dt 


( 12 ) 


= 0 


Momentum equation of gas 


, C \m(x,t) dp(x,t) n r, , n 
L Fr - + — + R Vl m(x,t) = 0 (13) 

dt dx 


where R rr is the regenerator flow resistance per unit length; 
L v , is the flow inductance per unit length; 

aev 1 

/? Fr = -—;L Fr = —, (14) 

/if r A fr 

where a = a + ( f3/d h )eR eh ; a = 175/(2eefh), j8 = 1 ,6/(2e 2 d h ) 
can be obtained from Tanaka’s results 5 ; d h = hydraulic 
diameter = edj( 1 - e); h = 033(k f /d h )R°* 7 based on 
Tanaka’s results 5 ; d m is the matrix wire diameter. 

For simplification, the axial variation of gas temperature 
and pressure at steady state are approximated by the aver¬ 
age values, i.e. T r (x) ~ T rm = (T u + T L )/2 and P,(x) ~ P tm 
= P ch . This was justified experimentally, as the temperature 
distribution in the regenerator is roughly linear. 

Equations (12) and (13) cannot be solved since the gas 
temperature T(x,t) is not known. The following energy 
equations for the gas and the screen matrix are thus derived 
using the above approximation on 7j.(x) and P r (x) 

Energy equation of gas 

C Fr dp(x,t) dm(x,t) 
y dt + dx 

+ Crr [T(x,t) - T s (x,t)\ = 0. (15) 

yXgr 


Energy equation of regenerator matrix 

T sr ^~^ + inx,t)-f(x,t)] = 0 (16) 


where C Fr and C Tr are the flow capacitance per unit length 
due to pressure change and temperature change, respect¬ 
ively; T Br and T sr are time constants of gas and matrix, 
respectively, and y = CJCy, x is the position measured from 
the hot side of the regenerator; 


n — ^ fr . /-■ _ ^rm^fr . 

^Fr n-r ’ Tr dt2 ’ 

i\i r m ”*rm 

_ P rm eV„C v _ p s (l - e)V 0 C s 

Tgr _ RT rm hA m ’ Tsr ~ hA m 


(17) 


where A HT is the surface area of the regenerator = 4V 0 ( 1 - 
e)/d m .d m is the wire diameter of screen disks; A' eh is the 
Reynolds number based on d h ; k f is the gas thermal conduc¬ 
tivity. 

Solution of Equations (12), (13), (15) and (16) can be 


obtained by Laplace transforms. Combining Equations (15) 
and (16) with (12) and (13), we obtain the gas continuity 
equation as 


dm(x,s) 

dx 


+ sCy r r p(X,S ) = 0 


and the gas momentum equation as 
dp(x.s) 


dx 


+ (R Fr + sLy)m(xx) = 0. 


(18) 


(19) 


where CpTr is the regenerator flow capacitance due to press¬ 
ure change and time responses of gas and matrix that is 
derived as 


CpTr — Cp 


1 + T s /[T gr ( 1 + s-rj] 
y + T s y[T gr ( 1 +ST sr )]' 


( 20 ) 


It is worth noting that Equations (18) and (19) are the 
wave equations. 

A linearly-perturbed dynamics model for the regenerator 
can be obtained 4 


Wr 0 (s) 

Pro(s) 


cosh( 


- - sinh(r r L r ) 

^cr 


- Z cr sinh( f r /, r ) cosh(r r L r ) 


x 


iKi(s) 
.Pd(s) - 


( 21 ) 


where r, = Vz r T r = regenerator propagation constant and 
Z cr = VZ,/T r = regenerator characteristic impedance; they 
are derived as 


F r — Ry r + .v/, Fr ); 

Zcr = (22) 

SC pT r 

Similar to the connecting tube, the dynamics model of 
the regenerator belongs to a distributed-parameter system. 
An equivalent circuit consisting of infinite sets of shunt and 
series impedances similar to Figure 3 can be drawn. This 
is based on the series expansion of cosh T, L, and sinh r,L, 
in Equation (21) with respect to r r L r . A system block dia¬ 
gram similar to Figure 4 can be used to represent the 
input/output relation of the regenerator. 

Equation (21) represents the dynamics model of the 
regenerator in terms of the pressure and mass flowrate that 
are derived from the conservation equations of continuity 
and momentum. The conservation equation of energy can 
also be simultaneously solved using the solutions of p r (x,s) 
and m T (x,s ) to obtain the gas temperature in the regenerator 
T fr (x,s). The gas temperature at the cold side of the regen¬ 
erator is finally derived as 


T f ,(L r ,s) = - 


Cpr T gr ( 1 + XT sr ) 
Ct r T, r 


Prof? ) 


y T gr ( 1 + .S'T sr )I j 

C Tr ST sr sinh[r r L r ] 

{m m (s) coshf T,i. r ] - niy,(.'•')}■ 


(23) 
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Figure 6 Equivalent circuit of displacer 


Displacer 


Linear network of Stirling refrigerator 


Taking force balance to the displacer will lead to a displacer 
dynamics model 

X d (s) = [A dw p„(s) - A de p c (s)] 

I _ 1 _ 

WjJ 2 + C„s + K, 

where K d is the spring constant, for gas spring at isothermal 
condition K d = p.AlgJV^: V gs = V ogs + A dgs X d C d is the 
coefficient resulting from the friction and the leakage 
losses; M d is the displacer mass. The equivalent circuit of 
the displacer is shown in Figure 6. 


(24) 


Expansion space 

The cold space is usually made of a small empty space, 
sometimes filled with a porous medium to enhance the heat 
transfer between the gas and the wall. The gas enters the 
cold space with phase difference between the temperature 
and the mass flow. The gas agitation in the cold space is 
so severe that an uniform temperature with T ce in the cold 
space can be assumed. Similar to the warm space, a lumped 
model is derived from the mass conservation equation 


m e {s) 


pjs) 

1 /sC e 


W*dc 

RT ce 


X d (s) 


(25) 


where C e = V C0 /(RT CC ) is the flow capacitance. The second 
term in the right-hand side of Equation (25) is the mass 
flowrate induced by the displacer motion, which is anal¬ 
ogous to a current source. The equivalent circuit of the 
expansion space is shown in Figure 7. 


rh e (s) 



Combining the equivalent circuits or block diagrams of all 
the components according to the process of the Stirling 
refrigerator, we obtain the linear network as shown in 
Figure 8. The block diagrams for the connecting tube 
and the regenerator in Figure 8 imply the equivalent cir¬ 
cuit of Figure 3. 

It is worth noting that there are two networks for the 
split-type free-displacer Stirling refrigerator: namely, fluid 
network and displacer network. The lower part of Figure 
8 is the fluid network with mass flow analogous to current 
and pressure analogous to voltage. The fluid network 
depicts the gaseous flow propagation from the compressor 
to the expansion space induced by the pressure. 

The upper part of Figure 8 is the displacer network 
which picks up the pressure difference signal (analogue to 
voltage) across the warm space and the expansion space as 
the input and the displacer motion X d (t) as the output. The 
displacer network thus depicts the displacer motion induced 
by the pressure force exerted on the displacer. 

Physically, the fluid network and the displacer network 
are coupled since the gaseous mass flow, the pressure, and 
the displacer have mutually affected each other. The sol¬ 
utions of gas flow, pressure and displacer motion can be 
determined analytically from the network analysis. 


Transfer function of split-type free-displacer 
Stirling refrigerator 

A block diagram can be drawn as in Figure 8 for the con¬ 
necting tube and regenerator with the distributed-parameter 
models. The transfer function model of the refrigerator can¬ 
not be derived directly from the circuit. Combining Equa¬ 
tions (1) to (24) and rearranging these will lead to an equi¬ 
valent impedance at the inlet of the connecting tube Z ti : 


Z ti (s) = 


G s (s) = R mm {s) 


~ E me (s) 


G g (j) + GMZqjs) tanh[r,(s)L,] 

G 7 (s) + fr;,f,v)/Z ct (,v)jtanh[I7.v)/,l 

W md (s)G 2 (s) 


1 + 


Rpm{s) + 


1 - G 5 (.s) 
G 2 (s)G 6 (s) 


1 - G 5 (s) 


(26) 


E md (s)G,(s) 

(27) 

1 -G 5 (s) 

G 7 (s) p to (s) = Gg(s)rh, 0 (s); 


G 0 (s) = R pm (s)W md (s) 

(28) 
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Figure 8 Linear network of split-type Stirling refrigerator 


G s (s) = G 2 (s)W md (sY, 

G 4 (s) = R pp (s) + /? pm (i)W m p(s) 

G 3 (s) = GM + G 2 (s)W Illp (*); 
G 2 (s) = D lic (s)R pm ( s ) 

G t (s) = D dw (s) + Z? de (s)R pp (s) 
sV„ 

E me (.v) - RT ~ 


W mp (s) = 
W md (5) = ' 


sV 

J v wo < 

*7V 

SP0A d» 
RT„ 


R mm (s) = R pp (s) = cosh(T r L r ); 
R pm G) = - Z cr sinh(T r T r ) 


D dc (s) = - 
O dw (s) = 


M d s 2 + C cl s + K d 

A ri 


dw 

M d s 2 + C d s + K, 


(29) 

(30) 

(31) 

(32) 


(33) 

(34) 


(35) 


Therefore, the linear network model in Figure 8 can be 
transformed to Figure 9, which is the equivalent circuit at 
the inlet of the connecting tube. 

The equivalent impedance at the current source m p is Z sys 
( 5 ) which follows 


1/Z sys (s) = 1/Z c (.v) + i/Z„ (s) (36) 

where Z c (s) = l/.vQ = RT</sV co . Therefore, the flowrate at 
the connecting tube inlet is 



Figure 9 Equivalent circuit at connecting tube inlet 


Substituting the above equations into the models of regen¬ 
erator and connecting tube and rearranging, we obtain a 
relation for 2f d (j) and X p (s) and another relation for p c (s) 
and X p (s). Two system transfer functions of the split-type 
free-displacer Stirling refrigerator are thus obtained 


Gdp(-s) 


^d(s) _ 

*p(i) 

' G : , (.v)G 8 (.v)/G 7 (.v) + G 2 (s) 
1 - G s (s) 


x 


Anp(^ ) Z n ( V) 4- f mm (T) 
1 + Z ti (s)/Z c (s) 


(38) 


m ti (j) = 


sp.AJiRFG - 
1 + Z„(.y)/Z c (.v) p 1 ’ 


(37) 


sPoA p 
RT C 
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G cp (s) 


= Pei*) 
x p (s) 


| G 4 (s) + 


G 3 (x)G 6 (x) 
1 - G 5 (s) 


Gg(j) 

G 7 (s) 


+ 


Rp mis) + 


G 2 (s)G b (s) 11 
1 - G,(s) JJ 


Net cooling capacity 

The net cooling capacity Q net can be evaluated by sub¬ 
tracting the heat losses from the maximum available coo- 
(39) ling capacity £> max . There are four types of heat loss: heat 
conduction loss of regenerator 0 cond , enthalpy flow loss of 
regenerator <2emh. shuttle heat loss of displacer <2 shllllle , and 
hysteresis loss of gas spring W ir , where 


fmp(T)^ti(s) fmmis) 

1 + Z ti (s)/Z c (s) 



£?shu — kf'nD.SlC-p 


T» - T l 
2cL r 


(44) 


where 

sinh( r,L.) 

Uj) = cosh(T t L t ); Lpis) = - -' (40) 

At 


where D e is the outside diameter of displacer; S,. is the dis¬ 
placer stroke; C x is the time fraction of heat transfer; c is 
the clearance between displacer and cylinder; L r is regener¬ 
ator length. Usually, C T = 0.2. The regenerator conduction 
loss is 


It is seen that the system input of the split-type free- 
displacer Stirling refrigerator is the piston displacement X p 
of the compressor. There are two system outputs: displacer 
displacement X d and expansion space pressure p c . The split- 
type Stirling refrigerator thus belongs to a single-input-mul¬ 
tiple-output (SIMO) dynamical system. 


Performance calculation 

Maximum available cooling capacity 

The maximum available cooling capacity of a free-displacer 
Stirling refrigerator can be evaluated by integrating the 
expansion space pressure p e (t) and volume V e (r), where 
V e (r) is related to the displacer displacement X d {t). For a 
given piston motion X p (t), X d (t) and p,X t ) can be obtained 
from the inverse of the system transfer functions G dp (s) and 
G cp (s), Equations (38) and (39). However, this compu¬ 
tation procedure is quite complicated. A simplification is 
needed. Since the piston is driven by a crank-shaft mech¬ 
anism, the piston and the displacer motion as well as the 
generated pressure and mass flow waves inside the refriger¬ 
ator approach sinusoidal. X d (t) and pjt) can thus be com¬ 
puted simply from the gain and the phase of the frequency 
response functions, G dp (ju>) and G ep (/'to). 

Assuming that a dp is the gain of G dp (/'w) with phase <f) d 
leading the piston; a ep is the gain of G cp (jw) with phase 
</> c leading the piston. Then 

Pc(t) =p pa sm(wt + </>e) 

= p eo sin(6>+ 4> e ) =p e (0) (41) 

X d (t) = X do sin ( wt + 4> d ) 

- X do sin( 0 + 4> d ) = X d (d) (42) 


Qcond = kef f ( ttDI/4) T " l T ' (45) 

k eff is the effective thermal conductivity of the regenerator. 
The enthalpy loss of the regenerator can be evaluated by 
integrating the mass flow rate and temperature waves at the 
boundary between the expansion space and the regenerator 

Genth =fj> mJt)C p T,( L,,t)dt 

= /7rC p |m ro | |7f(L r )| cos(<£ m - <K) (46) 

cf) m and are, respectively, the phases of m,J t) and 
Tf{L r ,t) leading the piston. The hysteresis loss of gas spring 
can be evaluated by the following relation 6 


W ir = \ y( y - 1)7’ w A dw ( A V og JV ogs ) 2 (47) 

The net cooling capacity of the refrigerator is calculated 
by the following relation: 

Gnet — Gmax — Gshu — Gcond ~ Genth — 14;: (48) 

The gas leakage loss at the clearance between the displacer 
and the cylinder wall is related to the seal design, manufac¬ 
turing process, and material used. It cannot be estimated 
accurately and is ignored. But the loss coefficient C d can 
be used to reflect this effect to some extent in addition to 
the displacer frictional loss. However, the gas leakage at 
the piston is ignored since it is small for a good com¬ 
pressor design. 


where 6 is the piston angle; p eo and X do are the amplitudes 
of p e (0 and XJt), respectively. The maximum available 
cooling capacity of the Stirling refrigerator is 

C2-n 

Gmax=/^ PedV, 

= /7ra ep a dp X 2 po A de sin( </> d - A ) (43) 

where X po is the amplitude of the piston, / is the 
operating frequency. 


System design analysis 

The system design analysis of a split-type free-displacer 
Stirling refrigerator can be carried out using the present 
linear network model. For simplification, the following 
assumptions are made in the analysis: 

1 Since the gas in the compressor and the connecting tube 
are assumed to undergo a stationary and isothermal pro¬ 
cess, the regenerator inlet temperature T r , approximates 
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displacer loss 
coefficient C d 




Figure 10 Flow chart of the linear network analysis 
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Pch=6.033 kg/cm~2 abs, f=35 Hz 
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Cd= 20 N s/in 
TL= 90 K 
Tc= 360 K 
Ttv= 330 K 

Displacer Stroke= 2.663 mm 
Displacer phase angle= 34.57 deg 
Expansion space pressure phase angle = -11.08 deg 
Gas Spring Constant= 2374 N/m, Wn/W= 1.341 
Displacer Natural Frequency fn= 46.93 Hz 
Net Cooling Capacjty= 1.143 \V 


2.5 


3.5 


4.5 


Expansion Space Volume, m~3 


5.5 

xlO -7 


Figure 11 PV diagram of the expansion space 
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Pch=6.033 kgf/cm~2 abs, f=35 Hz, CIST-1 Cooler 



Figure 12 Variation of Q net with f L 



Displacer natural frequency/operating frequency, fn/f 

Figure 13 Variation of Q„ e , with f 


the compression temperature T c and the connecting tube 
temperature T„ i.e. 

T n ~ T t ~ T c . (48) 

2 The structure of the warm space provides a high thermal 
inertia to the gas temperature. The gas temperature of 
the warm space was thus assumed constant and equal 
to the outside wall surface temperature T u and the 
regenerator inlet temperature 7’ rl . That is, 

T w » T ri « T h . (49) 

3 The gas temperature in the cold space T cc is approxi¬ 
mately equal to the outside wall surface temperature of 
the cold end, T t . This was verified experimentally 7 


since the convective heat transfer between the gas and 
the wall in the cold space is very large, especially for 
low T, and high operating frequency. That is, 

T ce - T l . (50) 

r L is assumed to be stationary since the thermal mass 
of the cold-end exchanger is large as compared to that 
of the gas inside the cold space. 

4 The mean pressure P is assumed to be equal to the 
charge pressure P ch . 

5 Since the tube wall usually has a large thermal mass 
with high thermal inertia, the wall temperatures at the 
two ends of the regenerator, T ri and T m , are approxi¬ 
mately constant or stationary and obey the following 
relations: 
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T, ~ 7' r , ~ T, « 7 C ; r ro « T ce ~ r L . (51 ) 


Given the refrigerator dimensions, the material physical 
properties, the operating conditions (l],P ch ,f,T c ,T^), the 
working fluid properties, and the value of C d , we can carry 
out system performance analysis according to the flow chart 
in Figure 10. The present design calculations are carried 
out for a test refrigerator with design specifications listed 
in Table 1. Shown in Figure 11 is the PV diagram of the 
expansion space. The variation of net cooling capacity Q net 
with the cold-end temperature T L is shown in Figure 12. 
<2 net is shown to increase with increasing T, and decreasing 
C d . The design calculation is also performed to investigate 
the variation of Q nel with the operating frequency /(= 
<x>/2tt). It is shown from Figure 13 that, for a fixed C d , 
Q nei first increases with increasing <o (or decreasing ot n /cu), 
reaches a peak value, and then decreases. An optimum Q„ a 
is obtained at a frequency ratio u>Jw around 1.2 to 1.5, 
depending on the value of C d . That is, the corresponding 
a> for the optimum Q nc[ ranges from 0.67 to 0.83w n . This 
coincides with the test results 8 - 9 . 

The frequency response functions of the refrigerator, 
G dp (/to) and G e p(/<n), can also be calculated. The results are 
shown in Figures 14 and 15. The amplitude and the phase 
of the displacer motion are shown to decrease with increas¬ 
ing operating frequency as shown in Figure 14. However, 
there is a lower limit (20 Hz for the simulated case) for the 
operating frequency below which <2 net will become negative 
(i.e. the refrigerator has no cooling capacity). For the gas 
pressure in the expansion space, the amplitude first 
increases with increasing frequency, reaches a maximum 
value then declines. An optimum amplitude exists (60 Hz 
in the simulated case). Though the maximum available coo¬ 
ling capacity Q max is proportional to the integral of p c and 
X d (i.e. the enclosed area in the PV diagram shown in Fig¬ 
ure 11), the operating frequency for the optimum p c may 
not result in an optimum Q„ n . This is caused by the vari¬ 
ation of the phase angle between p c and X d . 


Table 1 Design specifications of a test refrigerator 


1 . 

Compressor 



piston diameter 

26 mm 


piston stroke 

12 mm 


swept volume 

6.4 c.c. 

2. 

Connecting tube 



length 

190 mm 


inside diameter 

1.76 mm 

3. 

Warm space 



volume 

1 c.c. 

4. 

Gas spring chamber 



plunger diameter 

8 mm 

5. 

Regenerator 



wire screen material 

SS316 


wire screen mesh no. 

200 


wire screen length 

60 mm 


wire screen diameter 

11 mm 


number of screen disks 

880 

6. 

Displacer 



total mass 

27.3 g 


rod diameter 

12.7 mm 

7. 

Expansion space 



diameter 

12.7 mm 


Modification of linear network analysis 

Experimental verification of linear network 
analysis 

For verifying the linear network analysis, a test refrigerator 
was built in the present study. The design specifications are 
shown in Table 1. Helium gas with 99.999% purity was 
used as the working fluid. An electrical heating strip was 
taped on the surface of the cold head as the cooling load. 
The cold finger was put in a vacuum chamber for the 
measurement of net cooling capacity. The test results are 
shown in Figure 16. 

The deviation between the simulation results with C d = 
20 N m/s using the linear network analysis and the test 
results are small for T L < 105 K, but large (underestimated 
<2„e t ) for T, > 105 K. Figure 16 shows that the present 
linear network analysis with a constant C d is satisfactory at 
low 7 l . However, it will not be guaranteed all the time in 
selecting the value of C d , since C d depends on the displacer 
seal design, the manufacturing quality and the material 
used. There is no way for a designer to predict C d accu¬ 
rately without field experience. It may also not be true to 
assume a constant value of C d for various operating con¬ 
ditions. 

Modification of linear network analysis 

In order to improve accuracy in performance prediction, a 
modification of the linear network analysis is necessary. In 
the present analysis, a loss coefficient of the displacer C d 
is included. C d is the coefficient accounting for the friction 
and gas leakage losses of the displacer at the surface in 
contact with the cylinder of the cold finger and at the clear¬ 
ance seal of the gas spring. Hence, C d may not be constant 
but can vary with the displacer oscillating frequency / and 
cold-head temperature 7) . 

In the present study, the value of C d at various operating 
conditions is determined by fitting the test results of Q na 
with the calculated value. An empirical correlation of C d is 
obtained from the analysis of more than 200 test data, with 
T L in K and / in Hz 


1.0851 7.7346 x 10 6 

C d (/,T L ) = 430.66 + + 


(0.6653 + 0.2604/- 5.212 x 10” 4 / 2 ) (52) 


The analytical results using the above empirical relation 
coincide with the test results of Q„ el to within ±10%, as 
shown in Figure 16. C d obtained from the experiments 
includes the effect of gas leakage in addition to the fric¬ 
tional loss. 


Discussion and conclusion 

A linear network model is developed for the system per¬ 
formance analysis of split-type free-displacer Stirling 
refrigerators. The dynamics models are derived to describe 
the input/output relation of each component by use of the 
governing equations in conjunction with linearization and 
approximation. By connecting the equivalent circuits of the 
components together, a linear network consisting of a fluid 
network and a displacer network is obtained. Two transfer 
functions are derived for the displacer motion and the gas 
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Figure 14 Frequency response of G dp (s) 
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Figure 15 Frequency response of G ep (s) 

pressure of the expansion space. The system performance 
can then be evaluated from the frequency response of the 
two transfer functions by using the sinusoidal signal analy¬ 
sis. Since the analytical solutions are obtained, the compu¬ 
tation takes only a few seconds on a PC for each design cal¬ 
culation. 

Implementation of the linear network model in the sys¬ 
tem analysis of a split-type free-displacer Stirling refriger¬ 
ator relies on the correct selection of C d value. C d , however, 
is determined mainly from experience. The value of C d has 
been shown in the present study to vary with the operating 
frequency / and the cold-end temperature T h . An empirical 
correlation of C d is obtained. The performance prediction 
of a Stirling refrigerator using the empirical relation of C d 
is satisfactory. The C. d empirical relation can be further 


updated if more test results of various hardwares are avail¬ 
able. The present analysis is shown to be superior to the 
other methods in its high computation speed and good accu¬ 
racy. It is now possible for us to develop a powerful CAD 
tool for the design of split-type free-displacer Stirling 
refrigerators using the linear network model. 
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Pch= 6.033 kgf/cm-2 abs, f= 35 Hz 



cold—end temperature TL, K 

Figure 16 Comparison of the test results with the analysis 
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